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We] -known methods for solving the shape-from-shading problem require 



1^ ;d w len the reflectance map is not available, but is known to have a given 

^bme unknown parameters. This happens, for example, when the surface is 

La ni)ertian, but the direction to the light source is not known. We give 

algc :ithm that alternately estimates the surface shape and the light source 

Jse 



f the unit normal in parameterizing the reflectance map, rather than 

3tereographic coordinates, simplifies the analysis. Our approach also 

itera .ive scheme for computing shape from shading that adjusts the current 

the ocal normals toward or away from the direction of the Hght source. The 

idjui (ment is proportional to the current difference between the predicted 
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E, and a reflectance map, /?, the shape-from-shading problem may be 
of recovering a smooth surface, z{x,y), satisfying the image irradiance 



1 of E. Any given boundary conditions on z[x, y) should also be satisfied. 
es the form of a first-order partial differential equation (Horn, 1970 k, 



formulation are a number of assumptions, the principal one being that 

a surface patch does not depend on its position in space. Another is 

ipicts a smooth surface of homogeneous reflectance. Several algorithms 

to tackle the problem, notably those of Horn (1975), Strat (1979), and 
1981). 

difficulties these schemes face in practice is that the reflectance map 
lown. The reflectance map specifies how the brightness of a surface patch 
ientation under given circumstances. It therefore encodes information 
ing properties of the surface, and information about the distribution 
the light sources. In fact, the reflectance map can be computed from 
reflectance-distribution function and the light source arrangement, as 
; Sjoberg (1979). 

a new scene, the information required to determine the reflectance 
t available. Yet without this information, the sheipe- from- shad ing prob- 
mulatcd, much less solved. The dilemma may be resolved if a calibration 
shape appears in the scene, since the reflectance map can be computed 
lere we wish to consider the situation where we are not that fortunate, 
to evaluate how some basic assumptions can resolve this impasse, 
has looked at the problem of recovering shape from shading under the 
the image depicts a Lambertian surface illuminated by a point source 
5 unknown. Under the additional assumption that the surface is locally 
normals are shown to be recoverable by a local operation. This method 
on the iterative propagation of information across the image. 

serious drawbacks to the local approach, however. One problem is that 
ion is very restrictive. In fact, spheres are the only surfaces whose 
bifical. So this method naturally computes incorrect normals for other 
the errors for approximately spherical surfaces, such as ellipsoids of 
nay be acceptable. Further, the constraining effect of known occluding 
s cannot be incorporated into the local method. This is unfortunate 
mals provide powerful boundary conditions on the shape-from-shading 
n by Bruss (1983). Finally, because a local method does not take into 

when calculating the normal at a point, nearby normals may differ a 

arly in the presence of noise, 
it an alternative approach that does not sufifer from these disadvantages. 
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The task is to r icover the shape of a smooth surface depicted in an image, E, that is 
defined ove: a n gion H in the xy-plane. Let the shape of the surface be characterized 
by the fnnc :ion, n, that associates a unit normal with eeich point in Q. The problem 
is therefore to fi id n{x,y) over 0. Assume for now that the object has a Lambertian 
surface, anc tha it is illuminated by a single point source. If the vector s points to the 
source, and n(x, j) is the unit normal of a surface patch, then the ixpparent brightness of 
the patch is| give i by the reflectance map 

Rs{n{x,y)) =n{x,y) • s. 
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e way, force s to be a unit vector; this allows for the possibihty that the 
)urce nuiy be unknown. This way we can also deal with unknown sensor 
known surface albedo, provided it is uniform. 

then becomes one of finding a smooth shape, n, and source direction, s, 
mfige irradiance equation 

E{x,y) = n{x,y) • s \f{x,y) G Q. 
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tilce, b ightness caimot be determined with perfect accuracy, and so it appears 
nsform this into a minimization problem (Ikeuchi k Horn, 1981; Horn 
There is another reason to consider this as a minimization problem: 
to solve the image irradiance equation as it stands, we obtain a set of 
ions equivalent to the characteristic strip equations. Here, however, we 
;chcme lending itself to a parallel implementation on a grid, as originally 
n (1970). Further, a shape-from-shading problem that has noisy hnage 
t have a theoretical solution. A minimization approach will, however, 

rJjcovJry of a shape that fits the given data best, in a sense determined by the 



loth shape, n, and a source direction, s, that minimize 

/ / {E{x, y) - n[x, \j) • s) dx dy. 

s, and there are no errors in brightness measurements, then the image 
on will have been satisfied (although there is no guarantee that the 
c integrable; see Horn Sz Brooks, 1985). 
a regularizing component (Poggio & Torre, 1984) by incorporating the 
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to select a particularly smooth solution from a possibly infinite set of 
that a subscript here denotes partial differentiation, and that squaring 
cut to taking the dot-product with itself. Finally, we wish to insist that 
; length. This is accomplished with the constraint 



n'^{x,y) 



V(x,2/)en. 



2. All itrr<itiv( 



Coiiibining 



he 



which is to be 
relative imjiorta 
function iism to 
1985) 

Minimizing / 
and tliat / is tc 
the sohitiorils of 
functional 



has the Eul ;:r eq lation 



So, by subs 



where 



Now, a disc 



in which e is 
average of morniEis 



Rearranging 



complines s 



that 
Laplacian in 
ample, if we 



itiiti 



cIk 



t iree terms gives the composite functional 
f 



nimnnized with respect to n and s. Here, A is a scalar that weights the 

ce of the regularization term, while iJ,{x,ij) is a Lagrangian multiplier 

impose the constraint that n{x,y) be a unit vector (see Horn k Brooks, 

is a problem in the calculus of variations. First, assume that s is known 
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m, it follows that / has the Euler equation 
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pproximation to the Laplacian operator is given by 
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ape, given the light source direction. Other approximations for the 
to improved results, at the cost of increased computation. For ex- 
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/10> ). The simple 5-point approximation was adequate for our purposes. 

wc have yet to solve for /i, the hagrangian multiplier. This can he avoided, 
obs< rving that the division of the right hand side by (j f /i(c^/4A)) does not 
direc tion of the vector being computed. Since /i is intended to ensure that 
nor iialized, we can simply do this explicitly, as in 
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the problem of minimizing / with res])ect to s, given that n is known. 
1 in conventional calculus. Computing the partial derivative of / with 
lave 
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,he inverse of a jnairix. A discrete version of this formula, in which the 
iced by sums, is crisily obtained. An iterative scheme in both n and s 
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tage of the fact tliat s need not be a unit vector. If the briglitness of the 
we can normalize s so that it is a unit vector. Then the determination 
htly more complex, shice it involves a constrained minimization. 

ind performance of the scheme 



sc leme has two com])onents: one concerned with the recovery of shape, 
led with the determination of the source direction. The shape-recovery 
n intuitively satisfying form. In essence, a new normal is computed by 
erage, and adjusting this either tovv'ard or away from the source. The 
gn of the adjustment is determined by the brightness error of the current 
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lape, a new source direction is computed by a single pass through the 
pe-recovery, no iteration is necessary. The 3x3 matrix (nn'^') is summed 
as is the vector [En). The source direction can then be computed using 
tion or even Cramer's method (see Korn and Korn, 1908). The source- 
^nt has been tested on a number of images and sliapes. When the data 
the estimate of source direction is extremely accurate. Furthermore, 
very good in the face of significant noise. For example, a synthetic image 
a sphere illuminated by a point source in the direction (-4,3,8)'^'. The 
zed to 255 irradiance levels, and the correct surface normal was given 
250 image [)oints. Gaussian noise was added to the image giving an 
.ion in irradiance values of 34. Despite this, the source-finder computed 
urce direction that was only 2.7° in error. Further trials gave similar 



rection estimates are robust because the whole image is used. Theo- 
1cm is highly over-determined, as source direction is recoverable from 
of only three different surface orientations. Using the whole image 
that noise effects are significantly reduced. 

1-source-from sliading ])rol)lcm for a point light source and Lambertian 
ural two-way ambiguity. l( the image irradiance eciuation is satisfied 
ipe 111 «^nd the source direction Si, it will also be satisfied by the dual 
rce direction S2 where 
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viawmg direction. Here, both source direction and surface normals are 
t le viewing direction. This is easily verified by observing that 

np • ^|j =: 2z(ni • z) - ni) • (2z(si • z) - si) 

(z • z)(ni • z)(si • z) - 2(ni • z)(si • z) - 2(si • z)(ni • z) + rii • si 
1 si. 

init al values for the normals, the shape-and-source scheme will head for 
)thei of these solutions. The dual shape and source direction can then be 
imm idiately using the equations given above. 

pres iut two examples of the program at work. Each (synthetic) image con- 
ted a Lambertian surface illuminated by a point source in the direction 
e i] lages each contained more than 1000 points at which normals were to 
^ni|^|ed. ]|ormals were assigned an initial value of (0,0, 1)-^, as was the source di- 
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g boundary normals were given. The equations for n^y could be solved 



the Gauss-Seidel algorithm. Since we are ultimately interested in par- 

e^ilenta ion on a grid, we used the Jacobi method instead (despite the fact that 

3idel method has sHghtly better convergence properties). 

urn ge portrayed a hemisphere viewed from directly overhead. After 100 

th A = 0.005, the average angular difference between estimated and correct 

lesi than 3°. The maximum such deviation was less than 2.5 times the 

e. ' 'he estimate of the light source direction, at this time, had errors in 

zen th angles of 1.4° and 1.6° respectively. 

d n lage depicted a cylinder with rounded, hemispherical ends, viewed from 

erp( ndicular to its axis. After 60 iterations, this time with A = 0.003, the 

ar rror in surface normal was less than 5°. A further 30 iterations brought 

t > 4°. The maximum error remained somewhat larger, however, due to the 

to smooth the intersection between the cylinder emd the hemispheres. 

th and zenith angles for the source were 7.3° and 1.1° respectively, 

;r 90 iterations. These, too, improved slowly with further processing. 

t me \ as sometimes slow in converging. After rapid initial improvements, the 

] ess \ ould decrease appreciably. However, one might expect the scheme to be 

lome of the current iterative methods, given the disadvantage of not knowing 

irce direction. Convergence could be accelerated by employing multigrid 

t lat ] ropagate information across the image more quickly (see Terzopoulos, 

] estii gly, in the examples considered, a reasonable estimate for the source 

dfe obi ained after only a few iterations. Subsequent processing just improved 
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wo more new iterative schemes: the first extends the shape-and-source 
nations in which a simple model of the sky is also included; the second 
eloped above to find shape from shading, given a general reflectance 
recover source direction. 
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erative scheme for shape and source direction is composed of two parts. 
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